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ABSTRACT 

We have apphed two different automatic classifier algorithms to the BATSE Current 
GRB Catalog data and we obtain three different classes of GRBs. Our results confirm 
the existence of a third, intermediate class of GRBs, with mean duration ~ 25-50 s, as 
deduced from a cluster analysis and from a neural network algorithm. Our analyses im- 
ply longer durations than those found by Mukherjee et al. (1998) and Horvath (1998), 
whose intermediate class had durations ~ 2-10 s. From the neural network analysis no 
difference in hardness between the two longest classes is found, and from both meth- 
ods we find that the intermediate-duration class constitutes the most homogeneous 
sample of GRBs in its space distribution while the longest-duration class constitutes 
the most inhomogeneous one with <V/Vmax> ~ 0.1, being thus the deepest popula- 
tion of GRBs with Zmax ~ 10. The trend previously found in long bursts, of spatial 
inhomogeneity increasing with hardness, only holds for this new longest-duration class. 
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INTRODUCTION 

Since their discovery in the late 60's (Klebesadel, Strong & 
Olson 1973), GRBs have been a long-remaining puzzle (see 
Piran 2000; Meszaros 2001; Castro- Tirado 2001, for recent 
reviews). Models involving a short distance scale implied an 
emitted energy ~10''^ erg, whereas a cosmological origin re- 
quired 10^^ erg at least. The increasing degree of isotropy 
found as the sample of GRBs grew, plus the lack of faint 
GRBs (Meegan, Fishman & Wilson 1985) favored a cosmo- 
logical scenario. Moreover, after the launch, in 1991, of the 
BATSE instrument on board the CGRO, a very high de- 
gree of isotropy was found in the new, much larger sample. 
Finally, the measurement in 1997 of the first GRB redshift 
(Metzger et al. 1997) and the subsequent ones have con- 
firmed that GRBs are at cosmological distances. 

Concerning the physical mechanism of GRBs, the fire- 
ball shock model (Rees & Meszaros 1992; Meszaros & 
Rees 1993; Daigne & Mochkovitch 1998) is a progenitor- 
independent model for radiation emission that succeeds in 
explaining both the burst itself and its afterglow. There is, 
however, a variety of proposed objects that are capable of 
generating the GRBs (Nemiroff 1994) : from mergings of neu- 
tron stars with neutron stars (Paczyhski 1990) or with black 
holes (Narayan et al. 1992) to coUapsars (Woosley 1993; 
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MacFadyen & Woosley 1999) and hypernovae (Paczyhski 
1998). Magnetic instability in a neutron star being spun 
up by accretion in an X-ray binary could also produce 
them (Spruit 1999). There are more exotic models, involving 
quark stars (Ma & Xie 1996), mirror stars (Blinnikov 1999), 
or cosmic strings (Berezinsky, Hnatyk & Vilenkin 2001). In 
the face of such a boiling pot of theoretical ideas, to know 
the positions of the GRBs with respect to their host galax- 
ies and specially their redshifts are key issues (Bloom et 
al. 2000). Once a fair sample of redshifts becomes available, 
which should happen soon with the up-to-date technology of 
new missions like HETE II and Swift, the distances will be 
known and with that the released energy and the luminosity 
function of GRBs, together with their distribution across the 
Universe. That should certainly discriminate among existing 
models, and it should give as well unprecedented informa- 
tion on the very structure of the Universe up to redshifts far 
higher than 5 (Lamb & Reichart 2000; hereafter LROO), on 
the cosmic star formation history (Totani 1999), and on the 
evolution of galaxies (Totani 1997) . 

It is likely that more than strictly one progenitor could 
give rise to GRBs, since it has been shown that different ob- 
jects can produce a burst of gamma rays with the observed 
characteristics. Therefore, the catalog of GRBs may reflect 
the manifestations of various phenomena. To uncover them, 
various attempts to separate different classes of GRBs have 
been made. 

The two most generally accepted classes of GRBs are 
those arising from the bimodal distribution of their dura- 
tions (Kouveliotou et al. 1993; hereafter K93), that sepa- 



2 Balastegui et al. 



rates long (lasting for more than 2 s) and short (less than 
2 s) GRBs, the short bursts being at the same time spec- 
trally harder than the long bursts. The different spatial dis- 
tribution of the two classes has also been shown (Katz & 
Canel 1996) , and it is consistent with isotropy in both. How- 
ever, the longer GRBs appear to be more inhomogeneous in 
their space distribution than the shorter ones, as deduced 
from the higher value of <V/Vmax> for the short GRBs 
(this quantity measures the deviation of the space distri- 
bution from a homogeneous Euclidean distribution). Sepa- 
rating the long bursts class in two groups of hardness H32 
(the fluence ratio of spectral channel 3 to spectral chan- 
nel 2), respectively higher and lower than 3, Tavani (1998; 
hereafter T98) found that the long/hard bursts are more 
inhomogeneously distributed than the long/soft ones. We 
will show here that there exists, indeed, a trend, among the 
long bursts, of harder bursts being more inhomogeneous in 
their space distribution, which holds for hardnesses up to 
~4, with a slight rise of <V/Vmax> beyond that point. In 
the three-class classification that we propose here this trend, 
however, only exists for the longest class. 

In the following we first describe the two methods (clus- 
ter analysis and neural network) that are used to classify 
GRBs in the present work. The results from the two differ- 
ent methods are then discussed and compared. We find that 
both methods point to a classification in three classes, that 
resulting from splitting the 'old' long-burst class in two, and 
we compare the characteristics of the new three groups. Fi- 
nally, the possible physical meaning of the new classification 
is discussed. 

2 METHODOLOGY 

The usual approach in GRB classification has been based on 
the study of bivariate distributions. However, as noted by 
Bagoly et al. (1998; hereafter B98), the BATSE catalog pro- 
vides up to nine quantities intrinsic to the burst (7 related 
to energy and 2 related to duration), plus other quantities 
corresponding to spatial distribution and to errors in the 
magnitudes. New composite quantities can also be defined, 
such as the different measures of spectral hardness (from 
the fluence ratios in different spectral channels), and also 
V/Vmax. That involves a large number of variables which is 
difficult to handle, complex relationships among them (in- 
cluding nonlinear ones) being likely present. Such relation- 
ships can hardly be directly visualized, and thus multivariate 
analysis is needed. 

St arting from the BATSE Current GRB Catalog (avail- 
able at tittp : // www. batsc . msfc . nasa.gov/batse /grb/ catalog 
/current) in its version of September 2000, 1599 bursts have 
been selected: those for which nonzero values of all nine 
magnitudes are given. These magnitudes are the four time- 
integrated fluences Fchi-Fch4, respectively corresponding to 
the 20-50 keV, 50-100 keV, 100-300 keV, and 300+ keV spec- 
tral channels; the three peak fluxes P64, P256, and P1024, 
measured in 64, 256, and 1024 ms bins, respectively; and 
the two measures of burst duration T50 and T90, the times 
within which 50% and 90% of the flux arrives. Then a prin- 
cipal component analysis (PGA) of the standardized loga- 
rithms (zero mean and unity variance) of those quantities 
has been performed, obtaining results (Table 1) that are 
very similar to those of B98. As it is well known, PGA is a 



statistical method used in multivariate data analysis to ob- 
tain new variables, linear combinations of the original ones, 
which carry most of the variance of the system. Based on 
the correlations among the original variables, some of the 
new variables can be disregarded if they carry very little in- 
formation. For further details on PGA, the reader is referred 
to B98 and to Murtagh & Heck (1987; hereafter M87). 

As seen in Table 1, with only three variables, linear 
combinations of the original ones, ~96% of the system in- 
formation can be accounted for. The flrst row shows that 
64% of the information is carried by a variable which is a 
weighted sum of of all the original variables with nearly the 
same weight for each of them. The second principal compo- 
nent in importance is approximately the difference between 
the weighted sum of the logarithms of the three peak fluxes 
and that of the logarithms of the two durations, all again 
with similar weights. And with 5% of the total variance of 
the system the logarithm of the fluence in the fourth channel 
is found. 

Our current goal is to achieve an automatic classiflca- 
tion based on the nine original variables, and for that two 
different methods are used: a cluster analysis applied to the 
results of the PGA, and a neural network algorithm. 

2.1 Cluster Analysis 

For the cluster analysis the MIDAS statistical package has 
been used. As stated above, a PGA is flrst preformed. In 
this way are obtained new variables into which the problem 
becomes easier to separate. This result provides the starting 
point for the cluster analysis, where the Ward's criterion 
of minimum variance (Ward 1963; see also M87) is used. 
The analysis follows an agglomerative hierarchical cluster- 
ing procedure, which starts from n points spread over the 
9-dimensional space and groups them until ending up with a 
single cluster. The algorithm searches for clusters with min- 
imum variance among objects belonging to the same cluster 
and with maximum variance between clusters, and works 
with the 'center of gravity' of each cluster. That gives clus- 
ters as compact and as detached from each other as possible. 

A dendrogram is obtained, that shows the way groups 
are clustering, as well as the inner variance of the resulting 
groups. Thus, detecting a large rise in the variance by the 
union of two clusters means that two groups with remark- 
ably different characteristics have been merged. 

It is important to notice that the PGA looks for com- 
binations of variables to obtain axes which have the max- 
imum possible variance. As seen above, with just three of 
those axes more than 95% of the total variance is accounted 
for. The cluster analysis could, therefore, have been applied 
to those three new variables only. That has been checked 
and it is found that the results are barely affected by this 
reduction in the number of variables, since most of the sys- 
tem information is conserved. Here, however, the results of 
the complete 9-dimensional analysis will be presented, since 
those same nine variables will later be used for the neural 
network analysis. 

Goncerning the previous work of Mukherjee et al. (1998; 
hereafter M98), it must be noted that they choose six vari- 
ables for their analysis, three of them being the peak flux log 
P256 plus the two durations (log T50 and log T90), the other 
three being the total fluence log Ftotai and two hardnesses. 
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Table 1. Principal Component Analysis for the standardized logarithms of fluences, peak fluxes and durations. There is shown, 
in each row, the components of each principal axis in the base of our original variables (columns), together with the percentage 
of the variance carried by each of the new axes (first column). For instance, the first principal component is: — 0.39 log Fchl ~ 
0.40 log Fch2 - 0.40 log Fch3 - 0.33 log Fch4 - 0.22 log P64 - 0.28 log P256 - 0.36 log P1024 - 0.29 log T50 - 0.30 log T90. 



% 


log Fchi 


log Fch2 


log Fch3 


log Fch4 


log P64 


log P256 


log P1024 


log T50 


log T90 


64.3 


-0.39 


-0.40 


-0.40 


-0.33 


-0.22 


-0.28 


-0.36 


-0.29 


-0.30 


27.0 


+0.15 


+0.12 


+0.04 


-0.05 


-0.53 


-0.47 


-0.30 


+0.44 


+0.41 


4.9 


-0.22 


-0.19 


+0.06 


+0.92 


-0.10 


-0.13 


-0.17 


-0.06 


-0.08 


1.7 


+0.48 


+0.41 


+0.21 


+0.03 


-0.22 


-0.25 


-0.13 


-0.47 


-0.46 


0.8 


+0.56 


-0.05 


-0.77 


+0.19 


+0.16 


+0.06 


-0.09 


-0.01 


+0.12 


0.6 


+0.01 


+0.11 


+0.19 


-0.06 


+0.60 


-0.03 


-0.75 


+0.13 


-0.00 


0.4 


-0.02 


-0.05 


+0.11 


-0.02 


+0.04 


-0.04 


-0.07 


-0.69 


+0.71 


0.2 


+0.49 


-0.78 


+0.38 


-0.07 


-0.04 


+0.05 


-0.03 


+0.04 


-0.06 


0.1 


-0.01 


+0.08 


+0.01 


+0.00 


-0.46 


+0.78 


-0.40 


-0.03 


+0.01 



log H321 and log H32. It has been learned, from the PC A, 
that three variables are necessary, which together carry more 
than 95% of the variance, one of them being approximately 
the logarithm of the fluence in the fourth channel, Fch4 (see 
above). So it seems that M98 do lose information by not 
taking into account the fluence in the fourth channel sepa- 
rately, and also by not considering any peak flux later on. It 
must equally be noticed, when comparing results, that our 
GRB sample is twice as large as theirs. 

The main weakness of the cluster analysis is that it 
only deals with linear combinations of the variables. Such 
a weakness can be overcome by means of a neural network 
analysis, which also detects non-linear relationships. 



2.2 Neural Network 

Neural networks are artificial intelligence algorithms that 
can be used for an automatic and objective classification. We 
do not want to start from any prior classification. Therefore, 
a non-supervised algorithm is used. As we do not wish to in- 
troduce any tracer object either, the net is initialized at ran- 
dom. The 'Self-Organizing Map' algorithm (Kohonen 1990), 
implemented in the SOM_pak package from the Laboratory 
of Computer and Information Science of the University of 
Helsinki, is used. 

As in the cluster analysis, the entrance parameters are 
the logarithms of the same nine variables. 

The dimension of the output space must be specified, 
and based on the results of the cluster analysis the network 
is run two times asking first for a two and then for a three- 
dimensional output space, thus grouping either two or three 
classes of GRBs. The net is trained in two steps before look- 
ing for results. 



3 RESULTS 

3.1 Cluster Analysis 

In Fig. 1 the dendrogram with the last six levels of cluster- 
ing is shown. It can be seen that the first important increase 
of the variance occurs when joining group 3 with group 2, 
which tells that two groups with somewhat different charac- 
teristics have been merged, but the most significant rise in 
variance occurs when merging cluster 2 with cluster 1. From 
that it is concluded that there are two well-separated classes 
plus an emergent third class. 



Figure 1. Dendrogram of the 9-dimcnsional analysis. The num- 
bers at the bottom of the diagram arc identifiers of the groups, 
and those at the left show the inner variance of the groups. For 
instance, when merging group 6 with group 4 the variance of the 
cluster is 0.35. The largest increment in the variance occurs when 
merging group 2 with group 1, with a variance increase of 3.04. 

Fig. 2 shows what happens when adding, to the 
nine starting variables, the two extra variables H32 and 
<V/Vmax>. In that case the sample is reduced to 757 bursts 
only (instead of 1599), for which all the eleven quantities are 
known. It can be seen that the three-class classification is 
the most favoured one. 

Next, in Table 2, the main characteristics of each GRB 
class are shown, 2-1 and 2-II corresponding to the two- 
class classification, and 3-1, 3-II and 3-III corresponding to 
the three-class classification. The deviations correspond to 
a l\JN — 1. The results of the 11-dimensional cluster analy- 
sis are not shown here. They arc very similar to those of the 
9-dimensional one but are less significant because the sam- 
ple is reduced to one-half. We only comment that since the 
hardness has been added there, its weight has thus been en- 
forced, and then class I becomes slightly harder and shorter 
than when obtained from 9 variables. 

It must be noted that by just looking at the values of the 
dispersions in the four variables Tgo, H32, P1024 and Ftotai 
given in Table 2, it might seem that the variance would in- 
crease when shifting from the two-class to the three-class 
classification, but that is just an effect of projecting the 
groups onto these particular variables (three of them com- 
posites): the full 9-dimensional analysis shows the opposite, 
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Table 2. Characteristics of the classification from the 9-dimensional cluster analysis. T90 is in units of s, P1024 in photons cm ^s 

and Ftotal in units of 10^^' erg cm^^. 



Class 


N 


<T90> 


<H32> 


<V/Vmax> 


<Pl024> 


< Ftotal > 


<cos 9> 


<sin2 b-l/3> 


2-1 


580 


2.65±0.17 


5.96±0.20 


0.265±0.017 


1.29±0.08 


1.75±0.13 


-0.031±0.026 


-0.006±0.013 


2-II 


1019 


59.7±2.1 


3.11±0.05 


0.184±0.008 


3.33±0.20 


22.5±1.8 


-0.004±0.019 


+0.001±0.010 


3-1 


580 


2.65±0.17 


5.96±0.20 


0.265il:0.017 


1.29±0.08 


1.75±0.13 


-0.031±0.026 


-0.006±0.013 


3-II 


570 


51.3±2.3 


2.85±0.07 


0.296±0.012 


0.88±0.02 


4.58±0.21 


+0.021±0.025 


-0.000±0.013 


3-III 


449 


70.3±3.8 


3.43±0.06 


0.051±0.004 


6.44±0.41 


45.3±3.9 


-0.035±0.030 


+0.002±0.015 
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Figure 2. Dendrogram of the 11-dimcnsional analysis. Here 
the largest increase in the clusters variance occurs when joining 
groups 3 and 1, that suggesting three different classes of GRBs. 

as stated above. 

Adopting the same oxpoctod values for isotropy as for 
the 4B catalog, that is for the Galactic dipole moment <cos 
6> — —0.009 and for the quadrupole Galactic moment <sin^ 
b— 1/3> = —0.004, it can be seen in Table 2 that only one 
of the corresponding values for classes 3I-3III lies beyond 
1(7 of the expected value, and that is the dipole for the 3-II 
class, which is +1.2a above. Just such value being above 
la appears not significant and it is concluded that all three 
classes are isotropically distributed. 

In calculating the <V/Vmax> parameter, not all the 
1599 bursts could be used, but only those for which that 
value could be derived, and in a similar way when calcu- 
lating the dipole and quadrupole moments the GRBs that 
were overwrites (Fishman et al. 1999) were not taken into 
account. 

In Table 2 the class numbers are given in the order of in- 
creasing durations <T9o>. With the two-class classification 
the 'classical' GRB types are recovered: short/hard, which 
are fainter (taking as brightness the peak flux <Pio24>), 
and long/soft which are brighter and more non-Euclidean 
in their space distribution. As seen from Fig. 3, two classes 
with an overlapping distribution of durations have been ob- 
tained, in contrast with the classical definition of short, T90 
< 2 s, and long, Tijo > 2 s, GRBs. Now the short class has 
durations up to ^^20 s while the long-duration class starts at 
~2 s. This overlapping of the two classes was obviously sup- 
posed to exist, but based on the distribution of durations 
alone it could not be decided whether, in the overlapping 
region, a given GRB belonged to either of the two classes. 
Now the algorithm handles all the available magnitudes and 
assigns each GRB to the cluster to whose characteristics it is 



closer to. The hardness distribution (Fig. 4) does not difi^er 
significantly from that in K93. 

A first look at the three-class classification reveals that 
class I is exactly the same as in the two-class grouping: this 
is because the clustering method is agglomerative, which 
means that new groups axe formed by merging former ones, 
so the passage from three to two classes happens when merg- 
ing class II and class III GRBs. 

Let us concentrate on the new three-class classification. 
As stated already, 'old' class II has been divided into class 
II and class III. Class II is not properly an intermediate 
class: it has intermediate duration but still of the same order 
of magnitude as class III, and with an almost coincident 
distribution, as seen in Fig. 5. Class II is the softest and 
faintest class and the one most homogeneously distributed 
in space. Despite its duration being of the same order as 
that of class III, the fluence is one order of magnitude lower. 

The most striking result of this new classification is the 
extremely low value of <V/Vmax> in class III, which means 
that we are dealing with an extremely deep population that 
extends up to very high redshifts. Following the same proce- 
dure as in Mao & Paczyriski (1992), that is to calculate the 
theoretical value of <(Fniin/F)^''^> as a function of Zmax, 
and taking <V/Vmax> in Table 2 as an empirical value for 
that quantity, a value of Zmax for the distribution of the 
GRBs can be derived. Adopting a model Universe with Qm 
= 0.3, nA = 0.7, and Ho = 65 Km s'^ Mpc-\ the GRBs to 
be standard candles with a spectral slope a = 1 (Mallozzi et 
al. 1996), and assuming constant comoving GRB rate, one 
obtains: Zmax 
class II, and Zmax 
exceedingly high value, but £is will be seen in §4, class III 
can have very massive stars as progenitors, and in that case 
the GRB rate should be proportional to the star formation 
rate (SFR) rather than being constant. Taking as SFR(z) 
that of Madau & Pozzetti (2000), it is obtained for class III 
11.301^;^^, which is a more reasonable value. 
It was already known that separating long-class GRBs 
into two hardness groups results in two very different de- 
grees of inhomogeneity (T98). The <V/Vmax> values for 
the whole Current GRB Catalog have been calculated and 
the results are <V/Vmax> = 0.268±0.011, for GRBs with 
H32 < 3, and <V/Vmax> = 0.182±0.012, for GRBs with H32 
> 3. Now, with the three-class classification, very different 
degrees of inhomogeneity within the same interval of hard- 
ness are found, as it can be seen by comparing Table 2 with 
Fig. 6, and with little difference in mean hardness: those of 
class II and class III differ here with 4-10~^ significance in 
a Student test, but in the neural network classification the 
significance will only be of 0.42 while the <V/Vmax> dif- 
ference will remain. We leave a discussion of the value of 
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45.24+5;^'^ for class III. The latter is an 
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Figure 3. Duration distributions of classes 2-1 and 2-II. 
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Figure 4. Hardness distributions of classes 2-1 and 2-II. 
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Figure 5. Duration distributions of the three-class classification 
from cluster analysis. 

<V/Vmax> as related to hardness for §4. 

It can also be pointed out that in contrast with what 
happened in the binary long/short clEissification, where 
shorter bursts were harder, now class II is shorter than class 
III but it is slightly softer at the same time. 

3.2 Neural Network 

In the neural network case, how many classes arc to be ob- 
tained must be decided beforehand: knowing the dendro- 
grams that result from the cluster analysis, we ask for either 
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Figure 6. Hardness distributions of the three-class classification 
from cluster analysis. 



two or three classes. Their main characteristics are summar 
rized in Table 3, using the same units as in Table 2. 

As it should be expected, there are some differences in 
the composition of the classes as compared with those ob- 
tained from the clustering method, since the neural network 
method is not agglomerative. So, for instance, class I is no 
longer identical, in the two-group classification, to class I in 
the three-group scheme. Also, the 'short' GRBs which make 
up this class now have longer average durations than in the 
cluster analysis. 

There is also some change from the results of the cluster 
analysis in the three-group classification. Classes II and III 
now become more widely separated in duration, basically 
due to the decrease in duration of class II. The difference 
in hardness between class II and class III, in contrast, has 
decreased. 

As in the cluster analysis, all three classes are highly 
isotropic, with no value of the moments above O.Scr of the 
values expected for isotropy. 

Class II is now the intermediate class in peak flux. Prom 
both methods, cluster analysis and neural network, it is seen 
that despite the difference by one order of magnitude be- 
tween durations of class I and class II, their respective total 
fluences remain of the same order. 

The high inhomogeneity in the space distribution of 
class III is seen once more, and also how class II is again 
the most homogeneous one. Now the sample depths are: for 
constant comoving GRB rate, z^ax = 3.34lg;4|, 2.78t°]ll, 
and 15.451J 33 for classes I, II, and III respectively. For GRB 
rate proportional to the SFR, Zmax = 3.22^°' Jo, S.OSloJo, 
and G.GTly'gg- The value of Zmax ~ 11.3, obtained from the 



cluster analysis, corresponds to a Universe with an age of 
4.3-10* yr, and is in good agreement with the expectation of 
GRBs occuring out to at least z « 10. The value obtained 
from the neural network analysis of Zmax 6.7 corresponds 
to a Universe with an age of 8.8-10* yr. Both values for 
are below the redshift limit, z w 15-20, given for Population 
III stars by LROO. 

We conclude that the three classes respectively obtained 
from the cluster analysis and from the neural network algo- 
rithm show similar characteristics and thus both treatments 
are mutually consistent. 
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Table 3. Characteristics of the classification with the neural network. T90 is in units of s, P1024 in photons cm ^ s ^, and Fto^al 
in units of 10~^ erg crn^^. 



Class 


N 


<T90> 


<H32> 


<V/Vmax> 


<Pl024> 


<Ftotal> 


<cos d> 


<sin2 b-l/3> 


2-1 

2-II 


685 

914 


6.24±0.50 

63.5±2.3 


5.50±0.18 

3.12±0.05 


0.288±0.015 

0.159±0.008 


0.94±0.04 

3.82±0.22 


1.44±0.09 

25.1±2.0 


+0.002±0.024 

-0.024±0.021 


-0.005±0.012 

+0.001±0.010 


3-1 

3-II 

3-III 


531 
341 
727 


3.05±0.34 
25.0±1.4 
71.8±2.8 


6.20±0.22 
3.05±0.10 
3.15±0.05 


0.287±0.017 
0.307±0.019 
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Figure 7. Duration distributions of the three-class classification 
from neural network analysis. 
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Figure 8. Hardness distributions of the three-class classification 
from neural network analysis. 



4 DISCUSSION 

The two different automatic classifier methods above sug- 
gest the existence of three different groups of GRBs with 
different properties. It will now be examined whether such a 
classification docs make physical sense and we are actually 
dealing with three classes of GRBs. 

Let us first discuss the proposal by T98, of taking 
the long-burst class and dividing it into two groups with 
H32 higher and lower than 3 respectively. As seen in §3.1, 
long/hard bursts are more inhomogcneously distributed 
than long/soft bursts. This might seem to be in contradic- 
tion with the cosmological scenario, in which more distant 



bursts are expected to be softer due to the spectrum red- 
shift, and it leads us to conclude that class III bursts are 
intrinsically much harder than those of classes I and II. 

Fig. 9 shows that there is indeed evolution in hardness: 
the value of <V/Vmax> decreases with increasing H32. The 
hardness bins are taken so as to include similar numbers of 
bursts (~60) in each of them, in order to have comparable 
error bars. The value of <V/Vmax> is displayed in the po- 
sition of the mean of the hardness for each bin, and no error 
bars for the hardness are shown because the deviation is less 
than the symbol size, except for the last bin for which it is 
of about 0.2 

When <V/Vmax> decreases one is dealing with a more 
distant sample of GRBs, and then Fig. 9 tells us that, 
when sampling to higher distances, GRBs tend to be harder, 
and taking into account the hardness-intensity correlation 
(Dezalay et al. 1997) they should also be more luminous. 
This effect has to be interpreted, in a cosmological scenario, 
as a source evolution. There is a possible explanation: it is 
generally admitted that the upper limit of the stellar initial 
mass function (IMF) depends on mctallicity, and that lower 
metallicity allows more massive stars to form. When sam- 
pling GRBs farther away, one looks to a younger Universe, 
with lower metallicity, and thus with more massive stars. 
Therefore, if GRBs come from very massive stars those an- 
cient GRBs had sources with higher power and they were 
brighter and harder. 

Next are displayed, in Fig. 10 and Fig. 11, <V/Vmax> 
versus hardness for classes II and III, from the results of the 
cluster and of the neural network analyses, respectively. 

Here can be appreciated, with particular clarity in the 
Figure corresponding to the neural network result, how in 
the three-class scheme the 'new' class II no longer shows any 
trend of <V/Vmax> decreasing with increasing H32. Such a 
trend in the 'old' class II was due to the fusion into it of the 
'new' classes II and III, and now it is seen that class III is 
the only one to uphold the trend. To evaluate numerically 
the correlation, a Spearman rank test (see, for instance, its 
implementation in Press et al. 1992) has been applied, ob- 
taining for class III a Spearman-rank correlation coefBcient 
rs = —0.344, with a significance of 2-10~*, for the class re- 
sulting from the cluster analysis, and r^ = —0.354, with sig- 
nificance 4-10~^^, for that resulting from the neural network 
analysis. In contrast, for class II from the cluster analysis rg 
= 0.066 is obtained, with a significance of 0.26, and Ts = 
—0.051 with a 0.51 significance level for class II from the 
neural network. It is concluded, therefore, that class III re- 
ally shows clues of cosmological source evolution, which can 
be due to its being made of GRBs produced by very massive 
stars, likely through collapsars. 

Out of all bursts with known redshifts, eight of them 
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Figure 9. <V/Vmax> vs. Hardness for GRBs with T90 > 2 
s. The correlation between these two variables is clearly seen. 
Hardness bins have been taken so as to include similar numbers 
of bursts in all of them. Each hardness bin contains ~60 GRBs. 
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Figure 10. <V/Vmax> vs. Hardness for classes II and III from 
the cluster analysis. Here the nearly constant value of <V/Vmax> 
over the interval of hardnesses covered by class II is seen. The 
trend for class HI of lower <V/Vmax> with higher H32 can also 
be seen. As in Figure 9, the hardness bins have been taken so as 
to include similar numbers of bursts in all of them. 

entered into our classification. In the case of the neural net- 
work classification, seven out of those eight bursts are clas- 
sified as class III, and only one of them as class II. With 
the clustering method classification, six belong to class III 
and two to class II. One of these last two bursts is GRB 
980425, presumably related to SN 1998bw (Galama et al. 
1999), which is thought to be a peculiar GRB. In both clas- 
sifications GRB 970508 (Metzger et al. 1997) was assigned 
to class II. The assignation of any individual GRB to a given 
class by our methods is not entirely reliable, however, and 
has an uncertainty that is most important in the limiting 
region of each cluster. By looking at the scatter graphs, it 
has been checked that GRB 970508 lies in fact in the re- 
gion near class III. It must be stressed that, given the high 
isotropy of all three classes found here, there is no clue of 
Galactic structure for any of them. It could be axgued that, 
in the same way that no redshift has yet been measured for 
any GRB of class I, due to the fact that their being short 
makes them difficult to be detected with BeppoSAX (mainly 
sensitive to bursts longer than about 5-10 s), no redshifts of 
GRBs of class II have been measured either because they 
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Figure 11. <V/Vmax> vs. Hardness for classes II and HI from 
neural network analysis. In this graph, the trend for class III of 
lower <V/Vinax> with higher H32 is more evident than in Figure 
10. Same criterion as in Figure 10 for the hardness bins. 

are faint and their detection is equally hard. An alternative 
explanation would be that the fact that no afterglows from 
GRBs of classes I or II have been seen is rather due to their 
being produced by NS-NS or NS-BH mergings, which are 
expected to happen mostly outside galaxies, where the in- 
terstellar medium is too tenuous to produce any detectable 
afterglow (see Panaitescu, Kumar & Narayan 2001, for in- 
stance). That would also be consistent with the <V/Vmax> 
values found: while collapsars should appear first and be 
more frequent in the distant, early Universe, NS-NS and 
NS-BH mergings should start later and be more homoge- 
neously distributed down to low redshifts. In this context 
we can even speculate whether the differences between class 
I and class II GRBs might be due to one of them correspond- 
ing to BH-NS (or BH-WD) mergings, the other class being 
produced by NS-NS mergings. The whole question should 
be solved with the new generation of GRB detectors aboard 
the HETE II and Swift satellites. 

One should be concerned whether the structure of the 
GRB data may partially reflect instrumental biases. Hakkila 
et al. (2000) have suggested that the three-class classificar 
tion obtained by M98 might arise from a bias in measuring 
some bursts properties, such as duration and fluence, which 
would make some bursts in 'classical' class II to take the 
'new' class II characteristics (by lowering their duration and 
fluence). The fluence-duration bias, however, in spite of be- 
ing qualitatively understood, is not well quantified. Hakkila 
et al. based their analysis on M98 classes: their intermedi- 
ate duration class has durations Tgo of between 2 s and 10 
s, while the class II deduced here extends up to T90 longer 
than 100 s. Moreover, such bias acts on the farthest bursts, 
while what we find is that our 'new' class II GRBs are the 
closest ones. In addition, any bias that would make some 
bursts in the 'old' class II appear shorter and with lower flu- 
ence could hardly separate at the same time the evolutionary 
eflfects that we see in Fig. 9 into two groups: one with evo- 
lution (class HI) and the other one without it ('new' class 
II). 

5 CONCLUSION 

There are reasons to think that there exists more than one 
type of possible progenitors for GRBs, and each type may 
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give rise to groups of burst with different properties. We 
have searched for those groups in the current BATSE cata- 
log, with the aid of two automated classification algorithms, 
and confirmed that there exist two clearly separated classes 
of GRBs corresponding to the 'classical' classification of 
long/short GRBs. In addition, we have also obtained clear 
hints that there exists a third class, different from those pre- 
viously reported. An oversimplified way of looking at this 
would be to say that the third class arises from splitting the 
original long class into two groups with high and low peak 
fluxes respectively, in a similar way that the whole sample 
of GRBs has been divided, in previous studies, into pairs of 
groups according to duration (K93), hardness (T98), bright- 
ness (Nemiroff et al. 1994), or other characteristics (Pendle- 
ton et al. 1997). The present work, however, goes beyond 
that since nine quantities related to the bursts are used for 
the classification instead of taking a single parameter and 
then finding a value separating the bursts into two classes: 
there are overlapping zones in every original variable. What 
our procedures do is to trace a surface in the 9-dimensional 
space, separating classes from the way each variable relates 
to all others. Two different classes may well have the same 
duration or show nearly the same distribution for a given 
variable, but by taking into account the other variables as 
well these procedures still detect their existence. In contrast, 
univariate distributions would overlook them. 

Apart from the power of the method, the new group- 
ing of the bursts thus obtained has to be examined for its 
possible physical meaning and its correspondence with sepa- 
rate classes of GRB progenitors and/or mechanisms. Classes 
I, II, and III here defined correspond to different observa- 
tional depths (zmax) and may result from varying geometries 
of the observer with respect to the emitter, different param- 
eters of the explosion, or from different progenitors having 
different spatial distributions. Thus, every class has to be 
compared with several possible models. The physical sep- 
aration of classes II and III is strongly supported by the 
fact, which can hardly be due to chance alone, that having 
both classes together they show evolution of hardness and 
intensity with the maximum distance sampled while when 
separated such evolution only exists in class III. We con- 
clude, therefore, that class III, which likely has coUapsars as 
progenitors, is the one that can be detected up to very large 
redshifts, and it should thus be the most adequate one to 
learn about the history of the Universe at high z. We also 
suggest that classes I and II could correspond to NS-NS or 
NS-BH mergings instead. 
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